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ABSTRACT 


The  effects  of  thermal  nonequilibrium  (elevated  electron  tempera¬ 
ture)  and  ionizational  nonequilibrium  (finite -rate  recombination)  are 
studied  in  the  insulator  boundary  layer  of  a  potassium-seeded  nitrogen 
MHD  accelerator  flow.  The  nonsimilar,  compressible  boundary  layer 
is  assumed  steady,  laminar,  and  two-dimensional.  A  collisionless 
sheath  is  assumed  and  matched  with  the  boundary  layer  equations, 
which  are  solved  numerically  for  a  core  flow  Hall-neutralized  Faraday 
accelerator.  The  relative  importance  of  the  various  terms  in  the  elec¬ 
tron  energy  equation  are  assessed,  and  the  overall  effect  of  the  non¬ 
equilibrium  phenomena  on  the  boundary  layer  parameters  is  described. 
It  is  concluded  that  thermal  nonequilibrium  can  lead  to  significant 
B-wall  shorting  in  long  channels  and  that  Hall  effects  should  not  be  ne¬ 
glected  but  that  operation  is  not  noticeably  affected  by  ionizational  non¬ 
equilibrium  or  the  physics  of  the  electrostatic  sheath. 


iii 


A  EDC-TR -70-231 


CONTENTS 

Page 

ABSTRACT .  iii 

NOMENCLATURE .  vi 

I.  INTRODUCTION .  1 

II.  THEORETICAL  FORMULATION 

2. 1  Equations  of  Motion .  1 

2. 2  Gas  Properties .  5 

2.3  Boundary  Conditions .  10 

HI.  NUMERICAL  FORMULATION 

3.  1  Patankar-Spalding  Technique .  15 

3.  2  Pressure  Prediction  Technique .  19 

3.3  Initial  Profiles .  22 

IV.  RESULTS 

4.  1  Simultaneous  Thermal  and  Ionizational  Nonequilib- 

rium .  23 

4.  2  Electron  Thermal  Nonequilibrium .  25 

4.3  Ionizational  Nonequilibrium .  26 

4.4  Local  Hall -Neutralization .  26 

V.  CONCLUSIONS .  28 

REFERENCES .  29 

APPENDIX 

Illustrations 

Figure 

1.  Channel  Flow  Coordinate  System .  33 

2.  Control  Volume  for  Linearization  of  Conservation  Equa¬ 
tions  .  34 

3.  Non-MHD  Initial  Electron  Concentration  Profiles  for 

Different  Values  of  Recombination  Rate  Coefficient  at 

x/J20  =  0 .  35 

4.  Streamwise  Variation  of  Channel  Flow  Parameters  using 
Physically  Realistic  Recombination  Rate  and  Energy-Loss 
Factor 

a.  Geometry  and  Field  Configuration .  36 

b.  Core  Flow  Properties .  37 

c.  Boundary  Layer  Parameters .  38 

d.  Heat  Transfer  Rates .  39 


v 


AE  DC-T  R  -70-231 


Figure  Page 

5.  Development  of  Boundary  Layer  Profiles  using  Physi¬ 
cally  Realistic  Recombination  Rate  and  Energy-Loss 
Factor 

a.  Overall  and  Electron  Gas  Temperature  ....  40 


b.  Electron  Concentration  . .  41 

c.  Transverse  Current  Density .  42 

d.  Hall  Current  Density .  43 

e.  Velocity  . .  44 


6.  Comparison  of  Boundary  Layer  Profiles  at  x/iQ  =  82.  0, 
for  Different  Values  of  Energy-Loss  Factor 

a.  Electron  and  Overall  Gas  Temperature  ....  45 


b.  Electron  Concentration .  46 

c.  Transverse  Current .  47 

7.  Streamwise  Variation  of  Boundary  Layer  Parameters  for 
Different  Values  of  Energy-Loss  Factor 

a.  Displacement  Thickness .  48 

b.  Local  Heat  Flux .  49 

8.  Comparison  of  Boundary  Layer  Profiles  at  x/iQ  =  96.  5 
for  Different  Values  of  Recombination  Rate  Coefficient 

a.  Electron  Concentration .  50 

b.  Electron  Temperature .  51 

c.  Transverse  Current  Density .  52 


9.  Comparison  of  Boundary  Layer  Profiles  for  Core  Flow 
and  Local  Hall -Neutralization  at  x/iQ  =  35.2 

a.  Electron  and  Overall  Gas  Temperature  ....  53 


b.  Transverse  Current  Density  .  . .  54 

c.  Joule  Dissipation .  55 

10.  Streamwise  Variation  of  Boundary  Layer  Parameters  for 
Local  and  Core  Flow  Hall -Neutralization 

a.  Displacement  Thickness .  56 

b.  Local  Heat  Transfer .  57 


NOMENCLATURE 

A  Coefficient  in  linearized  equations,  given  in  Eq.  (60) 

A^2>,  A(5)  Coefficients  in  Eq.  (24),  defined  in  Ref.  10 
es  es 

a  Generalized  Conduction  coefficient  in  Eq.  (54) 
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Coefficient  in  linearized  equations,  given  by  Eq.  (61) 

Magnetic  field  intensity,  in  channel  flow  coordinate  system 
B  =  B(x)  k 

Generalized  source  term  in  Eq.  (54) 

Coefficient  in  linearized  equations,  given  by  Eq.  (62),  also 
used  for  generalized  conduction  coefficient  in  Eq.  (51) 

Saha  equilibrium  electron  concentration,  given  by  Eq.  (11) 

^  _  T 

Skin  friction  coefficient,  Cf  =  777; - 5- 

1  1/2  pug 

Potassium  neutral  concentration  at  Saha  equilibrium 

Original  potassium  seed  mass  concentration 
Specific  heat  of  overall  gas  at  constant  pressure 
Specific  heat  of  species  s  at  constant  pressure,  i.  e. , 
hs  =  CPS  Ts 

ps/p,  mass  concentration  of  species  s 

ws  -  v,  thermal  velocity  with  respect  to  mass  averaged 
velocity 

ws  -  vs,  the  instantaneous  velocity  of  a  particle  of  species  s 
with  respect  to  the  mean  velocity  of  that  species,  thus  the 
random  velocity  due  to  thermal  agitation,  i.  e. ,  ( <?g  ^  =  0 

dcidc2dc3,  an  element  of  velocity  space,  where  appropriate 
dw  or  dc  may  also  be  used  to  denote  an  element  of  velocity 
space 

dc^dc^dc^,  an  element  of  velocity  space 

dwidw2dw3,  an  element  of  velocity  space 

Electric  field  strength,  in  channel  flow  coordinate  system 
§  =  Exi  +  EyJ 

Absolute  charge  on  an  electron 

Non-normalized  velocity  distribution  function  for  s  species, 
i.  e. ,  ns  fs  dw 

Relative  speed  between  particle  r  and  particle  s 
(ai+l  a-j) 

(wi+l  "  ui-l)(ui+l  -  Ui) 
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(ai  +  aj.  i) 

(ui+1  -  «i_i)  (uj  -  Ui_l) 
h  +  v2/2,  stagnation  enthalpy 

L  Cshs>  static  enthalpy 
s 

(1/nig) (e^nt  +1/2  mgc|)+  p g/ps,  static  enthalpy  of 
species  s 

Unit  vectors  in  orthogonal  <x,  y,  coordinate  system, 

<~)  is  used  to  avoid  confusion  with  j,  the  current  density, 
and  with  i,  j,  and  k  when  used  as  dummy  indices 
— ► 

E  Js,  total  current  density 
s 

— ►  ^ 

Js  +  qsnsv,  total  current  due  to  species  s,  sum  of  con¬ 
ducted  current  and  convected  current 

eT,,  net  conducted  current 
s  s 

qsns\^s,  conduction  current  density  carried  by  species  2 

Equilibrium  constant,  njne/nK  =  2.41  x  10^  Z  e-50,408. /Te 

g  rs 

Boltzmann's  constant 
Electron  thermal  conductivity 

Channel  semi-width  in  direction  of  B-field,  a  function  of  x 
Reduced  mass 

Rest  mass  of  a  single  particle  of  species  s 

Number  density  of  particles  of  species  s,  i.  e. ,  ns  =  ffs  dw 

Prandtl  number 
(ui+1  -  Uj) 

4<xd  -  xxj)  (Uj+i  -  Ui-i) 

3 

4(xd  -  xu) 

(wj  -  Wj.i) 

4<xj)  -  xuHui+i  -  u>i-i) 

Scalar  pressure  p  =  E  p 

s  b 

3  x 

■g  ^  Ps{cici)  =  Jl  Ps(c2)s  =  nskTs'  pressure  of 

species  s,  which  is  defined  as  the  trace  of  |T_ 

s 
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u,  V,  w 


X,  y,  z 


Collision  cross  section  between  species  s  and  species  r 

E  q_,  total  heat  transfer 
s  ® 

ng  ^cs  (1/2  msc|  +  e?^.)}  heat  transfer  carried  by  species  s 

s/ 

Electric  charge  on  a  particle  of  species  s 

Overall  gas  constant 

Generalized  source  term  in  Eq.  (51) 

Temperature  of  the  heavier  components  of  the  gas 
Temperature  of  species  s 

Mass  averaged  velocity  in  x-,  y-,  and  z-directions, 
respectively 

vs  -  v,  average  drift  velocity  of  species  s 
Vg  =  vs  -?  =  v  +  <2)g  -  v  =  <c >s 


—  E  psvSJ  mass  averaged  velocity  v  =  ui  +  vj  +  wk  in  channel 
P  s 

coordinate  system 

Species  average  velocity,  i.  e. ,  vg  =  (w)g  =  ~  f?s  ws  dw 


Instantaneous  velocity  of  a  given  particle  with  respect  to 
laboratory  coordinates 


Channel  flow  cartesian  coordinate  system,  x  is  aligned  with 
the  primary  flow,  y  is  normal  to  the  electrode  walls,  and  z 
is  normal  to  the  insulator  walls 


Compressibility  factor 
Partition  function  =  Z^Ze/Zj^ 

Zrs  =  9.  11  x  10"1  +  7.  33  x  1CT5  T  -  1. 49  x  10-8  T2 


Recombination  rate  coefficient 


Hall  parameter,  j3  =  ce  B/nee 
-  nsVs  =  diffusion  flux  of  species  s 

Correction  factor  for  the  electrical  conductivity  of  a  multi - 
component  plasma,  which  takes  thermal  nonequilibrium  into 
account,  see  Eq.  (24) 

Potential  difference  across  sheath 
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Boundary  layer  thickness  at  u/uB  =  0.  95 
Effective  energy-loss  factor,  see  Eq.  (33) 
Displacement  boundary  layer  thickness 


A  multiplying  factor  included  in  the  energy  collision  integral 
to  account  for  inelastic  collisions;  for  perfectly  elastic 
particles  6g  =  1 

Permitivity  of  free  space 

Internal  energy  tied  up  in  an  individual  s  particle 

Recombination  energy  of  the  ions  of  species  s,  (per  particle)  - 
this  is  treated  in  the  formulation  of  the  equations  as  a  positive 
value 

Characteristic  vibration  temperature 
Defined  in  Eq.  (27) 

Debye  length,  see  Eq  (28) 

Electron-neutral  mean-free -path,  see  Eq.  (35) 

Overall  gas  viscosity 

Overall  electron  collision  frequency,  see  Eq.  (32) 


The  effective  momentum  collision  frequency  between 
particles  of  the  s  and  r  species, 

n  111m  nic;  n 

v *  = - —  ( - — - )/Q  c  3  f  dc' 

sr  3ngkTg  nip  +  ms  J  sr  s  s  s  . 

=  E  ps  =  mass  density 
s  ° 

=  msns  =  f  ms  fs  dw  =  mass  density  of  species  s 

(fys  J  ,  shorthand  notation  for  the  collision 

■®t  /chem 

term  in  the  species  conservation  equation 
Electrical  conductivity,  first  order  approximation 
Electrical  conductivity,  second  order  approximation 
Skin  friction 

Generalized  dependent  variable  in  Eq.  (51) 


x 


Stream  function 

Stream  function  at  centerplane  of  channel 
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0 
0j e 

01 

u 


a  parameter  in  Saha's  equation 

iVgk  1 

Normalized  stream  function,  0/0^,  see  Eq.  (52) 


SUBSCRIPTS  AND  SUPERSCRIPTS 


A,  B  Halfway  values  between  grid  line  i  and  the  neighboring 

grid  lines 

D,  U  Downstream  and  upstream  stations,  respectively 

i  i"1  component  or  a  grid  index 

o  Denotes  either  initial  value  at  channel  entrance,  or  value 

at  edge  of  sheath 

r  Denotes  some  species  other  than  s 


s  The  s  in  either  subscript  or  superscript  denotes  species, 

it  is  moved  up  to  the  superscript  position  when  there  are 
other  subscripts  (e,  electron;  i,  ion;  N2,  nitrogen) 

w  Denotes  wall  value 

x,  y,  or  z  Denotes  component  of  a  vector  in  the  x,  y,  or  z  direction 
m  Denotes  core  flow  value 
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SECTION  I 
INTRODUCTION 


In  recent  years,  a  great  deal  of  interest  has  been  focused  on  the 
effects  of  elevated  electron  temperature  and  nonequilibrium  recombi¬ 
nation  on  the  performance  of  magnetohydrodynamic  (MKD)  accelera¬ 
tors  and  generators.  The  problem  is  particularly  acute  for  accelera¬ 
tors,  where  higher  current  densities  tend  to  elevate  the  electron  tem¬ 
perature,  which  in  turn  drastically  affects  the  electrical  conductivity. 
Hale  and  Kerrebrock  (Ref.  1)  and  Sherman  and  Reshotko  (Ref.  2)  have 
examined  the  effects  of  elevated  electron  temperatures  on  insulator 
walls  assuming  equilibrium  electron  concentration  throughout  the 
boundary  layer,  as  determined  by  Saha's  equation  evaluated  at  the  elec¬ 
tron  temperature.  For  practical  accelerators  and  generators,  this 
Saha  equilibrium  assumption  is  adequate  in  the  core  flow  and  across 
most  of  the  boundary  layer.  Near  the  wall,  however,  there  is  evidence 
which  suggests  that  the  electron  concentration  is  governed  by  finite  re- 
combintation  rates.  In  this  case,  the  very  strong  coupling  between 
electron  temperature  and  electrical  conductivity,  through  Saha's  equa¬ 
tion,  is  considerably  modified.  If  this  ionizational  nonequilibrium  ex¬ 
tends  very  far  out  into  the  boundary  |layer,  then  the  effect  on  the  overall 
boundary  layer  might  be  significantly  different  than  earlier  analyses 
indicate. 

The  purpose  of  this  study  is  to  provide  insight  into  the  relative  im¬ 
portance  of  elevated  electron  temperature  and  nonequilibrium  recombi¬ 
nation  in  typical  accelerators.  The  flow  in  a  typical  Faraday  accelera¬ 
tor  is  studied,  using  potassium  seeded  nitrogen  as  a  working  fluid. 

Only  the  insulator,  or  B-wall,  is  considered,  since  this  is  where  the 
most  potentially  detrimental  effects  are  expected  to  occur.  To  make 
the  problem  more  tractable,  steady,  two-dimensional,  laminar  flow  is 
assumed  in  the  plane  of  the  magnetic  field.  A  nonsimilar,  implicit, 
numerical  technique  is  used  to  calculate  the  flow  field  in  the  continuum 
region,  and  this  is  carefully  coupled  to  a  relatively  simple  collisionless 
sheath  model. 


SECTION  II 

THEORETICAL  FORMULATION 
2.1  EQUATIONS  OF  MOTION 

The  procedure  followed  in  the  derivation  of  the  equations  of  motion 
was  the  conventional  multifluid  formulation.  The  details  of  this  approach 
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are  presented  in  several  standard  references,  but  the  most  nearly 
applicable  formulation,  for  the  present  purpose,  is  contained  in  Sutton 
and  Sherman  (Ref.  3).  This  formulation  briefly  consists  of  taking  the 
mass,  momentum,  and  energy  moments  of  the  Boltzmann  equation  for 
each  species,  thus  yielding  macroscopic  equations  describing  the  species 
behavior.  This  very  general  set  of  multifluid  equations  is  then  simplified 
in  Ref.  4  to  a  two-fluid  system  of  equations  of  boundary  layer  form.  It  is 
believed  that  this  two-fluid  system  of  equations  adequately  represents  the 
nonequilibrium  phenomena  under  study,  while  offering  the  advantage  of 
being  simple  enough  to  be  numerically  tractable. 

The  major  simplifying  assumptions  inherent  in  this  analysis  are 
summarized  as  follows: 

1.  Magnetic  Reynolds  number  small  (Maxwell's  equations 
ignored), 

2.  Plasma  quasi-neutral  (except  in  sheath), 

3.  All  species  except  electrons  in  thermal  equilibrium, 

4.  Usual  steady,  laminar,  two-dimensional  boundary  layer 
approximations  valid, 

5.  Wall  cold,  fully  catalytic,  and  a  perfect  electrical  insula¬ 
tor, 

6.  Debye  length  smaller  than  or  equal  to  the  electron  mean- 
free -path, 

7.  Three-dimensional  swirl  due  to  Hall  current  neglected, 

8.  Three -body  recombination  only,  and 

9.  Radiation  effects  negligible. 

The  two-fluid  formulation  of  the  equations  of  motion  is  presented 
as  follows: 


Overall  Continuity 

9(pu)  9  ( pw) 

-  +  -  =  0 

9x  9z 

Overall  Momentum 

9u  9u  dp  9  /  9u\ 

pu  —  +pw  —  - - +  j  B  +  —  (u  —  ) 

9x  9z,  dx  *  9z  \  9z/ 


(1) 


(2) 
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Overall  Energy 


3H 

3H 

3 

P 

3H 

—  + 

pw  —  = 

— 

3x 

3z 

3z 

Pr 

3z 

3 

5k 

E 

+  j  E 

— 

(T 

X  X 

Jy  y 

3z 

.2 

3z 


l\  y  3  (u  ) 


Pr/  2  3z 


^  1 

rec 


3z 


3T 

k  —Z 
e  3z 


Electron  Diffusion 


3C  3C 

pu  — —  +  pw  — — 


3x 


3z 


3 

3z 


Ve. 


+  P. 


Electron  Energy 
3T 


3T  13 

e  e 

pu  — -  +  pw  — —  - - 


3x 


3z  C  3z 
e 


k  n  3T 
_ e_ _ e 

cp  32 

l  *e 


m  re  3T 
e  ez  _ e 

C  3z 
e 


,  ExJex  +  EyJey  _  l  f^\  ^ 


CeCP, 


5  \m  /  a 
'  e  f  o 


6eff  <Te  ”  T) 


2e‘ 


T  + 
e 


rec 


5k 


Pe  +  ep  u^x 
ce  raeCPeao 


(3) 


(4) 


(5) 


Equations  (1)  through  (5)  describe  the  boundary  layer  on  the  insu¬ 
lator  wall  of  a  typical  accelerator  channel,  as  illustrated  in  Fig.  1. 
The  primary  flow  is  in  the  x-direction,  the  imposed  magnetic  field 
(B-field)  is  in  the  z -direction,  normal  to  the  insulator  wall,  the  pri¬ 
mary  electrical  current  flow  is  in  the  y-direction,  and  the  Hall 
current  flows  in  the  x-direction.  Equations  (1)  through  (3)  are  basi¬ 
cally  the  conventional  boundary  layer  equations.  In  this  two-fluid 
formulation,  they  describe  the  motion  of  the  heavy  components  of  the 
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gas,  namely  the  nitrogen  background  gas,  the  potassium  neutrals,  and 
the  potassium  ions.  In  addition  to  the  conventional  boundary  layer 
terms,  they  contain  terms  accounting  for  the  MHD  body  force,  the  joule 
dissipation,  the  energy  carried  by  electrons  and  ions  diffusing  down 
through  the  boundary  layer  toward  the  cold  wall,  and  the  energy  con¬ 
ducted  by  the  electron  gas  due  to  the  electron  temperature  gradient. 

Equations  (4)  and  (5)  are  boundary  layer  equations  which  describe 
the  electron  gas  boundary  layer.  Equation  (4),  the  electron  diffusion 
equation,  is  familiar  to  those  who  have  worked  with  high  temperature 
boundary  layers  involving  species  diffusion.  The  electron  energy 
Eq.  (5),  by  contrast,  is  not  in  such  common  use.  In  the  mathematical 
formulation  of  Hale  and  Kerrebrock  (Ref.  1),  for  example,  the  joule 
dissipation  term  was  simply  balanced  off  against  the  collision  term. 

This  was  quite  adequate  in  the  free  stream  and  outer  boundary  layer 
but  may  be  inadequate  near  the  cold  wall  where  the  conduction  and  dif¬ 
fusion  terms  dominate.  The  electron  energy  equation  used  by  Sherman 
and  Reshotko  (Ref.  2)  was  analogous  to  Eq.  (5)  in  degree  of  sophistica¬ 
tion,  but  this  sophistication  was  to  some  extent  negated  by  their  assump¬ 
tion  of  ionizational  equilibrium. 

Equations  (1)  through  (5)  include  terms  such  as  the  current  densi¬ 
ties  jx  and  jy,  the  electron  diffusion  flux  rez,  the  electron  recombina¬ 
tion  rate  pe,  the  pressure  p,  and  the  equilibrium  electron  concentration 
Cg.  These  are  related  to  the  dependent  variables  by  the  following  aux¬ 
iliary  relations: 


Ohm's  Law 
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R  e  combination 


P„  =  a 


rec 
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(9) 


State 


p  =  pZRT 
Saha  Equilibrium 
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The  notation  used  in  the  preceding  equations  is  defined  in  the  Nomen¬ 
clature,  while  expressions  are  presented  in  Section  2.  2  for  the  gas 
properties  M  ,  Pr,  Z,  R,  Cp,  tce,  arec,  6eff*  CTe*  and  04. 

Equations  (6)  and  (7)  are  Ohm's  Law  with  the  assumption  that  Ex 
is  a  constant  through  the  boundary  layer  and  has  a  value  such  that 
jx  =  0  in  the  free  stream;  i.  e.  ,  the  Hall  current  is  neutralized  in  the 
free  stream.  Equation  (8)  is  the  ambipolar  diffusion  flux  in  the  form 
given  in  Ref.  4.  Equations  (9),  (10),  and  (11)  can  be  found,  for  exam¬ 
ple,  in  Ref.  3. 

The  foregoing  system  of  five  partial  differential  equations,  together 
with  the  auxiliary  relations,  the  gas  property  relations  in  Section  2.  2, 
and  the  one -dimensional  continuity  equation,  forms  a  closed  system  in 
terms  of  the  dependent  variables,  u,  w,  H,  Ce,  and  Te,  and  the  inde¬ 
pendent  variables  x  and  z. 


2.2  GAS  PROPERTIES 

Accurate  evaluation  of  the  gas  properties  is  one  of  the  most  impor¬ 
tant  prerequisites  to  successful  MHD  flow  calculations.  This  is  because 
all  of  the  gas  properties  are  affected  by  the  current  density  j.  The  most 
strongly  affected  is  the  electrical  conductivity  ae,  while  j  itself  is  pro¬ 
portional  to  oe.  Thus  any  inaccuracy  in  the  calculation  of  cre  will  feed 
back  through  j,  to  multiply  the  inaccuracy.  For  this  reason  a  great  deal 
of  effort  has  been  made  to  utilize  the  best  available  methods  to  calculate 
the  gas  properties. 
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The  overall  gas  properties  R,  Pr,  Z,  Cp  and  n  are  not  greatly 
affected  by  the  potassium  seed  and,  therefore,  are  essentially  the 
same  as  for  pure  nitrogen.  R  and  Pr  were  assumed  constant  at  values 
of  1776  ft3/sec3-°R  and  0.  75,  respectively.  The  following  expressions 
were  used  for  the  remaining  overall  gas  properties: 

*3.661  x  10"3  T*| . 4964  -  -55-'--31- 

Z  =  1  +  276.6  -  e  T  (12) 

P 


r  ,  ibm  i  r  t  -|°-668 

y  =  1.21  x  10~'>  -  -  (14) 

ft-secj  L300°Kj 

where  T  is  in  degrees  Kelvin  and  p  is  in  atmospheres.  Equations  (12) 
and  (13)  are  due  to  Ring  of  ARO,  Inc. ,  and  are  valid  within  a  pressure 
range  from  0.1  to  1.0  atmosphere  and  a  temperature  range  from  200  to 
5000°K.  Equation  (14)  is  a  curve  fit  to  calculations  by  Baulknight  (Ref.  5) 
and  is  accurate  to  within  3  percent  from  300  up  to  5000°K. 

The  three-body  recombination  reaction 

K  +  e""  K+  +  2  e”  (15) 


is  considered  to  be  the  dominant  reaction  in  this  analysis.  The  kinetics 
of  this  reaction  have  been  studied  by  Curry  (Ref.  6),  Cool  and  Zukoski 
(Ref.  7),.  and  Dugan  (Ref.  8)  for  potassium  and  other  alkali  metals. 
Calculations  from  these  references  have  been  distilled  by  Demetriades 
of  STD  Corporation  into  the  curve  fit 


=  3.47  x  10"20  T  -4‘765 
rec  e 
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for  the  recombination  rate  coefficient  c^gg,  where  orec  is  in  rn^/sec 
and  Te  is  in  degrees  Kelvin.  This  curve  fit  agrees  well  with  Refs.  6, 
7,  and  8  up  to  3500°K  and  closely  approximates  the  results  of  Ref.  8 
up  to  8000°K. 

The  conductivity  of  a  gas  for  species  s  is  defined  in  th-  notation  of 
Sutton  and  Sherman  (Ref.  3)  as 


where  mgr  is  the  reduced  mass 
collision  frequency. 

By  evaluating  Eq.  (17)  first  for  the  potassium  ions,  the  dominant 
resistance  to  the  diffusion  of  ions  through  the  partially  ionized  plasma 
is  assumed  to  be  collisions  with  diatomic  nitrogen  molecules.  Accord¬ 
ingly,  the  summation  in  Eq.  (17)  is  approximated  as 

E  n>*r  v*r  .  m*^  nN2  <c,iN^>  <QiN.,)  (18) 
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ns 


l  m 


(17) 
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sr  sr 
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and  i/sr  is  the  effective  momentum 


For  Maxwellian  velocity  distributions  of  the  ions  and  neutrals  about 
the  heavy  species  temperature,  the  average  relative  speed  (g)  between 
colliding  particles  is 


(19) 


Weber  and  Tempelmeyer  (Ref.  9)  present  integrated  potassium  ion- 
nitrogen  cross  sections  which  may  be  approximated  by  the  curve  fit 


~  2.39 


x  l<f17  t'1-756  x  10'5  T 


0.5) 


where  T  is  in  degrees  Kelvin  and  QiN9  is  in  meters^.  Equations 
(17)  through  (20)  then  combine  to  yield  the  ion  conductivity 


(20) 


ffi  = 


2.39  x  10 
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Evaluating  Eq.  (17)  for  the  electrons  is  not  nearly  so  simple  be¬ 
cause,  for  partially  ionized  plasmas,  no  single  interaction  predomi¬ 
nates.  Thus  it  is  necessary  to  account  for  electron -neutral  collisions, 
electron-ion  collisions,  and  electron-electron  interactions.  This  has 
been  done  in  a  very  general  formulation  by  Demetriades  and  Argyropoulos 
(Ref.  10)  over  the  entire  range  of  molar  ionization  fraction  above  10“?. 
This  expression  is 
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(22) 
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The  foregoing  expressions  are  only  as  good  as  the  cross  sections, 
Qeg,  which  are  used  in  them.  In  this  study,  the  integrated  collision 
cross  sections  used  by  Garrison  (Ref.  11)  are  adopted  as  follows: 

QeK  -  4.0  x  io"18  (29) 

Q_M  =  4.5565  x  10’20  +  2.81786  x  10"23T  -  4.99704 x  1027T2 
“  ™  ^  ©  6 


+  3.30643  x  10-31T3  (30) 

e  i 

where  Qe£  and  QeN3  31:46  111  meters^,  and  Te  is  in  degrees  Kelvin. 
Equation  (30)  is  a  curve  fit  based  on  data  collected  by  Shkarofsky, 
Bachynski,  and  Johnston  (Ref.  12).  Garrison  (Ref.  11)  has  experi¬ 
mentally  verified  the  accuracy  of  Eq.  (22)  when  used  with  Eqs.  (29) 
and  (30)  for  near -equilibrium  potassium-seeded  nitrogen. 


The  electron  thermal  conductivity  Ke  can  be  expressed  from  the 
results  of  Demetriades  and  Argyropoulos  (Ref.  10)  as 


(31) 


where 
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8k  T  \ 1/2 

^r) 


—  n  Q  . 
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(32) 


and  cre,  arQ,  Qei,  Q6k>  and  Q6n2  are  given  in  Eqs.  (22),  (23),  (26),  (29), 
and  (30),  respectively. 


The  average  energy-loss  factor,  which  accounts  for  inelastic 

collisions  in  the  collision  term  of  Eq.  (5),  may  be  defined  in  the  notation 
of  Sutton  and  Sherman  (Ref.  3)  as 
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where  6g  is  the  energy -loss  factor  for  the  s  species,  as  contrasted  to 
6eff  which  is  averaged  over  all  the  species.  For  perfect  elastic  colli¬ 
sions  between  electrons  and  particles  of  species  s,  6S  =  1.0.  The  use 
of  5eff  offers  the  advantage  that  it  can  be  evaluated  experimentally  with¬ 
out  detailed  knowledge  of  the  various  species  densities  and  the  momen¬ 
tum  cross  sections.  Experiments  have  been  carried  out  by  Demetriades 
(Ref.  13),  with  diatomic  nitrogen,  which  indicate  that  6eff  remains 
approximately  constant  at  about  7.  0  x  10"^,  for  N2  temperatures  from 
about  2800  to  5000°K.  This  corresponds  to  a  alaout  18.  While 

the  conditions  of  this  experiment  were  by  no  means  identical  with  those 
of  the  boundary  layer  problem  under  consideration,  it  is  nevertheless 
felt  that  this  is  as  reliable  as  any  data  currently  available. 

It  should  be  pointed  out  that  all  of  the  electron  gas  properties  are 
strongly  dependent  on  empirical  evaluation  of  the  collision  cross  sections 
or  energy -loss  factors;  thus  any  future  work  undertaken  should  be  based 
on  improved  evaluations  of  these  quantities  as  they  become  available. 


2.3  BOUNDARY  CONDITIONS 

In  channel  entrance  flows  such  as  the  one  considered  here,  the  free- 
stream  boundary  conditions  on  u,  H,  Ce,  and  Te  can  be  calculated  by 
simply  dropping  the  transverse  terms  from  the  linearized  equations  of 
motion.  In  this  way,  the  downstream  values  of  the  dependent  variables 
can  be  calculated  explicitly  in  a  way  that  is  consistent  with  the  boundary 
layer  calculation,  which  is  implicit  in  the  z -direction  and  explicit  in  the 
x- direction. 

The  wall  boundary  conditions  are  straightforward  for  the  overall 
gas  and  are  specified  as 

!u  0 

w  -*■  0  (34) 

H-  Hw(x) 

The  wall  boundary  conditions  for  the  electron  gas  properties,  Ce  and  Te, 
cannot  be  specified  in  a  straightforward  way  because  of  the  presence  of 
the  electrostatic  sheath.  This  is  considered  in  the  next  section. 
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2.3.1  Sheath  Analysis 


The  physics  of  the  electrostatic  sheath  is  determined  by  the  follow¬ 
ing  assumptions: 

1.  Wall  fully  catalytic,  perfect  electrical  insulator, 

2.  Debye  length  <  electron  mean-free-path,  (collisionless 
sheath),  and 

3.  Electron  and  ion  distributions  Maxwellian  about  their 
respective  temperatures. 


Of  the  foregoing  sheath  assumptions,  the  second  is  the  most  restrictive 
for  the  present  purposes.  The  electron  mean-free-path  Ae  near  the 
wall  may  be  approximated  by 


(35) 


while  the  Debye  length  A^  is  given  by  Eq.  (28).  These  two  equations 
may  be  combined  to  yield 


Ad 
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where  A(j/Ae  <  1  from  assumption  number  two.  In  a  typical  potassium- 
seeded  nitrogen  accelerator  with  inlet  core  flow  conditions  of 

p  s  1 .  atm 

T  =  Te  ^  3500°K  (37) 


the  prevailing  conditions  near  the  insulator  wall  might  reasonably  be 
expected  to  be  about 
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(38) 


Tg  s  700. °K 
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Evaluating  Eq.  (36)  at  the  conditions  of  Eq.  (38),  using  the  cross  section 
given  by  Equation (30),  yields  X(j/Xe«  1.  Thus,  the  collisionless  sheath 
assumption  is  somewhat  marginal,  and  if  subsequent  calculations  yield 
electron  concentrations  significantly  smaller  than  lO"^,  this  assump¬ 
tion  will  be  invalidated. 


The  fully  catalytic  wall  assumption  means  physically  that  no  elec¬ 
trons  or  ions  that  reach  the  wall  will  bounce  back  without  recombining. 
Since  the  wall  is  a  perfect  electrical  insulator,  the  sheath  potential  will 
adjust  itself  such  that  equal  fluxes  of  electrons  and  ions  reach  the  wall. 
The  thickness  of  the  sheath  is  roughly  one  Debye  length,  the  electron 
mean-free -path  is  assumed  greater,  hence  the  sheath  is  relatively 
collisionless. 


The  electron  flux  to  the  wall  may  be  expressed  as 


re 


z 


dc^dc2dc2 


where  the  limits  on  c^  and  C2  are  -  ®  to  + 
to  -t2eA0/mg]^^.  Since  fg  is  Maxwellian, 
ion 


=  -  n 
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®,  but  the  limits  on  C3  are  -  ® 
this  integrates  to  the  express - 

e  A  cf) 

“  ]r  <V " 
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The  ion  flux  integral  is  similar,  except  that  the  limits  on  C3  are  -  ® 
to  0,  yielding  for  the  ion  flux 


ri  =  “  n0 
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k  T 


2im-i 


1/2 


(40) 


Since  the  sheath  potential  drop,  A<£  ,  seeks  a  value  which  will  equate 
the  electron  and  ion  flux,  it  may  be  expressed  as 


A  4 


k  T 
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(41) 


Combining  Eqs.  (39)  and  (41)  yields  the  following  expression  for  the 
electron  flux  across  the  sheath: 
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where  the  subscript  o  dentotes  that  a  quantity  is  evaluated  at  the  edge 
of  the  sheath.  For  the  heavy  gas  properties  these  values  are  practi¬ 
cally  equal  to  the  wall  values,  but  for  the  electron  gas  they  may  be 
many  orders  of  magnitude  different  from  the  wall  values. 


Equation  (42)  is  one  of  two  expressions  which  are  used  to  describe 
the  effect  of  the  sheath.  The  other  expression  is  the  heat  flux  which  is 
carried  across  the  sheath  by  the  electron  gas.  It  is  obtained  in  a  simi¬ 
lar  way  by  expressing  the  electron  energy  flux  as 


=  j  (I  meCe) 


fe  dc1dc2dc3 


with  the  same  limits  as  before.  This  integrates  to 


= - ZZ 
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Inclusion  of  Eq.  (41)  yields  the  final  form  of  the  electron  heat  flux 
across  the  sheath 


Pp  Cep  k  Teo 


k  T  1 1/ 2 
_ o 

2  wik 


(44) 


Equations  (42)  and  (44)  may  be  used  to  define  boundary  conditions 
on  Ce  and  Te,  if  they  .are  handled  in  such  a  way  as  to  be  compatible 
with  the  continuum  boundary  layer  equations.  This  is  accomplished 
through  the  use  of  a  transition  region,  the  physics  of  which  is  dis¬ 
cussed  in  the  next  section. 

2.3.2  Sheath-Continuum  Compatibility 

The  transition  region  is  herein  defined  as  the  region  bordered  by 
the  collisionless  electrostatic  sheath  on  the  wall  side  and  the  collision- 
dominated  continuum  flow  which  extends  throughout  most  of  the  bound¬ 
ary  layer  and  all  of  the  core  flow.  The  transition  region  is  physically 
thin,  extending  from  a  normal  distance  of  about  one  Debye  length  from 
the  wall,  to  a  distance  of  perhaps  a  thousand  electron  mean-free-paths. 
Each  species  within  the  transition  region  is  assumed  to  be  Maxwellian 
about  its  own  temperature,  ccilisional  energy  transfer  between  species 
becomes  important  with  increasing  normal  distance,  but  electron- 
electron-ion  collisions  are  infrequent  enough  to  render  recombination 
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effects  negligible.  The  neutral  gas  properties  are  assumed  constant 
at  their  wall  values  across  the  transition  region,  and  convection  terms 
are  neglected.  Current  densities  are  assumed  low  enough  for  joule 
dissipation  to  be  neglected. 


With  these  assumptions,  Eqs.  (4)  and  (5),  the  electron  diffusion 
and  energy  equations,  become 


—  k]  - 

dz  L 


(45) 


and 


d3e. 


dz 


3  C, 


ep 


2  C 


-  Seff  (Te  -  T> 


(46) 


where  re is  given  by  the  ambipolar  diffusion  expression,  Eq.  (8),  and 

dT^ 

(47) 


qe  =  -  *  -  +  nt  h  rp 

*  e  dz  e  e  % 


When  the  boundary  conditions  Ce  =  CeQ  and  -  reQ  at  z  =  0  are 
given,  Eq.  (45)  can  be  integrated  assuming  that  ore  »  cr^  and 

dT  dC 
— ®  ® 

dz  dz 


This  yields  the  following  expression  for  the  Ce  profile  across  the  tran¬ 
sition  region: 

C  «  Ce  +  m  z  (48) 

e  o 

where 

*2  c. 

n  =  -  re  - - -  (49) 

o  k  (T  +  T)cr. 

c  X 

which  is  approximately  constant  since  Ce  from  Eq.  (21),  and  Te 
does  not  vary  markedly  from  TeQ  across  the  transition  region. 

Equation  (46),  which  describes  the  electron  temperature  variation 
across  the  transition  region,  is  somewhat  more  difficult  to  handle,  but 
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given  the  boundary  conditions  Te  =  TeQ  and  qGz  =  qeo  at  z  =  0,  it  can  be 
formally  integrated  as  follows: 


This  expression  is  integrated  numerically  across  the  transition  region 
to  yield  values  of  Te  and  dTe/dz  at  the  interface  between  the  transition 
region  and  the  collision-dominated  outer  boundary  layer. 


By  using  assumed  values  for  CeQ  and  Tgo,  Eqs.  (48)  and  (50)  are 

used  to  evaluate  Ce,  Te,  and  their  derivatives  where  the  transition 
region  merges  into  the  boundary  layer  flow.  These  values  then  serve 
as  mixed  wall  boundary  conditions  for  Eqs.  (4)  and  (5),  the  electron 
gas  boundary  layer  equations. 


The  numerical  method  of  solving  the  equations  of  motion  coupled 
with  the  sheath  dominated  transition  region  equations  is  described  in 
the  next  section. 


SECTION  III 

NUMERICAL  FORMULATION 


3.1  PATANKAR-SPALDING  TECHNIQUE 

Equations  (1)  through  (5),  the  boundary  layer  equations,  were 
solved  simultaneously  using  a  numerical  technique  developed  originally 
for  turbulent  flow  by  Patankar  and  Spalding  (Refs.  14  and  15).  This 
technique  was  chosen  for  several  reasons,  chief  of  which  was  the  fact 
that  it  would  facilitate  the  later  extension  of  the  solutions  into  the 
turbulent  flow  regime.  Only  the  essential  features  of  the  Patankar- 
Spalding  technique  are  discussed  here,  with  the  modifications  that  are 
required  for  this  particular  application.  A  much  more  general  treat¬ 
ment  is  to  be  found  in  Ref.  14,  whereas  Ref.  15  is  recommended  for 
the  details  of  the  computer  program. 
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The  Patankar -Spalding  technique  takes  advantage  of  the  fact  that 
Eqs.  (1)  through  (5)  share  the  common  boundary  layer  form 


3$  3  d  3 

Pu  —  +  pw  —  =  — 
3x  32  3z 


+  S 


(51) 


where  <f>  is  a  generalized  dependent  variable,  u,  H,  Ce,  or  Te.  This 
system  of  equations  is  solved  in  the  normalized  von  Mises  coordinate 
system,  where  u  is  the  cross -stream  coordinate,  defined  by 


dio 


pw  pu 

—  dx  +  —  dz 


h 


(52) 


where  ^  is  the  center  plane  stream  function,  a  constant,  defined  by 


cl 

d*  =  J  pu  dz 
0 


(53) 


In  the  transformed  (x-u)  plane  Eq.  (1),  the  continuity  equation,  is  auto 
matically  satisfied,  and  Eq.  (51)  assumes  the  form 


3(|> 

3x 


+  b 


(54) 


2 

where  a  =  pn/ip^C  and  b=  S/pu.  The  expressions  for  a  and  b  are  given 

in  Appendix  D  of  Ref.  4  for  the  boundary  layer  equations,  Eqs.  (2) 
through  (5). 


By  following  Patankar  and  Spalding  (Ref.  15),  Eq.  (54)  is  linearized 
by  considering  the  shaded  control  volume  in  Fig.  2;  d  is  assumed  to  vary 
linearly  in  u  between  the  u  grid  lines,  with  all  d  values  inside  the  control 
volume  evaluated  at  the  downstream  station.  With  these  assumptions  in 
mind,  the  convection  term  is  approximated  as 


3d 


rxD  rwB  3d 
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where  Pj,  P2,  and  P3  are  defined  in  the  Nomenclature.  The  conduc- 


tion  term  can  be 

approximated  as 

aD  - 

_U  +i+l 

AD  -  *D 
_U  *i-l 

3 

3<j> 

aB 

ui+l  "  “i 

aA 

“i  "  Wi-1 

3ui 

cl 

3w 

1  . 

z  g5  (^i+1  *i)  "  g6  (*i  ^i-l)  (56) 


where  g5  and  gs  are  defined  in  the  Nomenclature.  The  source  term, 
b,  is  in  most  cases  approximated  by  the  expression 


(57) 


Under  some  circumstances,  it  is  desirable  to  evaluate  a  particularly 
important  source  term  to  a  higher  degree  of  accuracy.  Patankar  and 
Spalding  (Ref.  15)  have  suggested  that  this  be  done  by  extending  the 
evaluation  of  the  term  to  the  adjacent  two  w  grid  points.  Here  the 
source  term  b^  is  assumed  to  vary  linearly  in  w  between  adjacent 
u  grid  lines.  The  averaged  b°  is  then  evaluated  similarly  to  the  con¬ 
vection  term  as  follows: 


du  dx 


^XD  “  XU^  ^UB 


w. 


) 


Plbi+1  +  P2bi  +  P3bi-1 


<*D  '  V 


(58) 


Equation  (57)  can  then  be  substituted  into  Eq.  (58)  for  an  improved 
approximation  to  the  source  term, 

It  should  be  pointed  out  at  this  time  that  there  are  no  hard  and  fast 
rules  for  the  evaluation  of  source  terms.  Equations  (57)  and  (58)  are 
excellent  guidelines,  but  their  application  often  requires  considerable 
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ingenuity.  For  example,  the  ionization  term'  &e  in  the  electron  diffu¬ 
sion  Eq.  (4)  varies  several  orders  of  magnitude  between  adjacent 
w-grid  lines  within  the  boundary  layer.  Since  it  is  the  dominant  term 
in  the  electron  diffusion  equation  throughout  most  of  the  boundary  layer, 
one  might  consider  using  Eq.  (58),  supposedly  a  finer  approximation 
than  Eq.  (57),  to  evaluate  it.  As  it  turns  out,  Eq.  (58)  yields  a  con¬ 
siderably  poorer  approximation  due  t=-  the  pact  that  it  is  dominated  by 

the  b^?+1  term.  Thus  for  exponentially  varying  source  terms,  a  step¬ 
wise  approximation  is  better  than  a  xx.iear  one. 

The  linearized  version  of  Eq.  (54)  may  be  obtained  by  substituting 
the  approximate  expressions  Eqs.  (55),  (56),  and  (57)  into  Eq.  (54), 
which  yields,  after  some  manipulation, 

♦i  *  A*?+l  +  B+?-l  +  c  (59) 


where 


g*  -  P- 


A  - 
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(62) 


where  the  coefficients  A,  B,  and  C  are  evaluated  at  the  upstream 
station.  Equation  (58)  could  also  have  been  used,  resulting  in  slightly 
different  coefficients  in  Eq.  (59). 

Once  the  boundary  layer  equations  are  linearized  to  the  form  of 
Eq.  (59),  it  may  be  observed  that  the  coefficients  form  a  tridiagonal 
matrix.  These  matrices  may  be  inverted  in  a  straightforward  way  by 
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successive  substitution  to  yield  each  succeeding  downstream  profile 
in  <j> .  This  procedure  is  described  in  detail  in  Patankar  and  Spalding 
(Ref.  15). 

One  of  the  most  attractive  features  of  the  Patankar -Spalding  tech¬ 
nique  is  its  use  of  Couette-flow  assumptions  near  the  wall.  As 
Patankar  and  Spalding  conceived  it,  the  Couette-flow  region  is  a  com¬ 
putational  necessity  when  used  with  their  turbulent  flow  hypotheses. 

For  laminar  flow  of  the  overall  gas,  the  Couette-flow  model  is  not 
necessary,  as  very  little  accuracy  is  lost  by  solving  the  linearized 
equations  all  the  way  to  the  wall.  For  the  electron  gas  boundary 
layer,  however,  the  Couette-flow  concept  offers  a  method  of  accurately 
coupling  the  mixed  sheath  boundary  conditions  to  the  linearized  electron 
diffusion  and  energy  equations.  Here  the  transition  region  of  Section  2.  3 
is  treated  as  a  Couette-flow  region,  and  an  iteration  is  performed 
between  the  Couette-flow  solutions,  Eqs.  (48)  and  (50),  and  Ce  and 
Te  boundary  layer  profiles.  The  iteration  proceeds  as  follows:  first 
CgQ  and  TeQ  are  set  equal  to  their  upstream  values.  Then  Eqs.  (48) 
and  (50)  are  solved  to  yield  Te,  dTe/dz,  Ce,  and  dCe/dz  at  the  outer 
edge  of  the  transition,  or  Couette-flow,  region.  In  general,  these  quan¬ 
tities  are  not  equal  to  Te,  dTe/dz,  Ce,  and  dCe/dz,  evaluated  at  U2.5' 
and  the  slip  values  and  wall  values  of  Ce  and  Te  are  readjusted  to  bring 
them  into  equality.  This  readjustment  yields  new  values  for  Ceo  and 
T0o,  and  the  iteration  is  repeated  until  C6o  and  Teo  no  longer  change 

in  successive  iterations. 

While  this  iteration  is  cumbersome  and  expensive  in  computer  time, 
it  offers  the  distinct  advantage  that  the  boundary  conditions  need  not  be 
matched  at  the  edge  of  the  sheath,  where  continuum  flow  is  a  poor 
assumption.  Rather,  the  boundary  conditions  are  matched  at  something 
like  a  thousand  electron  mean-free -paths  out  from  the  wall,  where  the 
Patankar-Spalding  technique  is  accurate,  and  the  flow  is  still  essentially 
continuum  one -dimensional  Couette-flow. 


3.2  PRESSURE  PREDICTION  TECHNIQUE 

In  order  to  calculate  the  downstream  velocity  profiles,  a  pressure 
gradient,  dp/dx,  must  be  specified  for  the  overall  momentum  equation, 
Eq.  (2).  Since  the  flow  is  confined,  the  pressure  gradient  cannot  be 
specified  a  priori  for  a  given  area  variation.  This  difficulty  can  be 
surmounted  by  several  types  of  numerical  schemes,  all  of  which  involve 
either  assumptions  regarding  the  variation  of  dp/dx,  or  iteration  between 
dp/dx  and  the  overall  equations  of  motion.  The  latter  techniques  were 
rejected  in  this  study  because  of  the  excessive  computation  time  they 
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involve.  Two  methods  of  the  former  type  were  used,  the  first  for  a 
priori  specified  channel  area  variation,  and  the  second  for  a  priori 
specified  temperature  variation. 

3.2.1  Specified  Area  Variation 

A  method  due  to  Kitowski  (Ref.  16)  utilizes  the  upstream  pres¬ 
sure  gradient  dp/  dx  in  the  overall  momentum  equation.  Once  the  down¬ 
stream  velocity  u  and  temperature  T  profiles  are  calculated,  the  down¬ 
stream  pressure  is  calculated  using  the  integral  continuity  equation  as 
follows: 

pu  p  u 

dw  ■  —  dz  = - dz  (63) 

R  ZT 

where  the  equation  of  state,  Eq.  (10),  and  the  definition  of  the  norma¬ 
lized  stream  function,  Eq.  (52),  have  been  included.  Equation  (63)  can 
be  integrated  from  the  wall  to  the  center  plane  as  follows: 

R  fl  ZT 

p  =  -  I  —  dm  (64) 

A (x)  J0  u 

to  yield  the  downstream  pressure  as  an  integral  function  at  the  down¬ 
stream  u  and  T  versus  to  profiles,  and  the  specified  downstream  channel 
semi-width  f(x).  Once  this  downstream  pressure  is  known,  Eq.  (63) 
may  again  be  integrated  to  yield  the  z  value  corresponding  to  each  u 
value;  i.  e. , 


z  = 


(65) 


Kitowski' s  method,  although  it  fails  to  rigorously  satisfy  the  momen¬ 
tum  equation,  yields  satisfactory  results  for  the  flows  under  consideration 
here. 


3.2.2  Specified  Temperature  Variation 

Often  it  is  desirable  to  design  a  channel  configuration  based  on  the 
a  priori  specification  of  certain  properties  down  the  channel.  The  case 
of  near  constant  core  flow  temperature  and  current  density  is  considered 
in  detail  in  this  investigation,  since  typical  accelerator  designs  often 
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follow  this  process.  For  constant  static  temperature,  Eqs.  (2)  and  (3), 
the  overall  momentum  and  energy  equations,  may  be  reduced  in  the  core 
flow  to 

du  dp 

pu  —  =  -  —  +  j  B  (66) 

dx  dx  Y 

and 

,  du 

pu  —  =  j  E  +  j  E  (67) 

x  x  y  y 


Equations  (66)  and  (67)  may  then  be  combined  to  yield 
dp 

u  —  =  -  j  E  +  j  (uB  -  E  ) 
dx  xx  Y  y 
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therefore, 

dp  +  j2 

_  _  Jx  Jy 

dx  u  o 

e 


(68) 


Equation  (68)  specifies  the  pressure  gradient  which  is  required  for  con¬ 
stant  temperature  one -dimensional  flow.  However,  various  anomalies 
arise  in  practice  because  of  the  presence  of  the  rather  thick  boundary 
layers.  These  anomalies  are  accounted  for  by  simply  multiplying 
Eq.  (68)  by  a  correction  factor  (T/TQ)  to  an  arbitrarily  large  power, 
thus 

.  Ox  +  jv) 

dx  u  a 

e 

where  TQ  is  the  static  temperature  at  the  channel  inlet.  Then  as  the 
computation  proceeds  down  the  channel,  if  the  temperature  drifts  away 
from  Tq  the  pressure  gradient  is  modified  to  keep  the  temperature 
nearly  constant. 
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With  the  pressure  gradient  specified  by  Eq.  (69),  the  downstream 
pressure  is 


PD 


<XD  "  V 


(70) 


which  yields  sufficient  information  to  evaluate  Eq.  (65)  for  the  z  values 
corresponding  to  each  w.  It  is  useful  to  use  this  method  of  specified 
temperature  variation  to  design  a  channel,  and  then  proceed  with  the 
method  of  specified  area  variation  to  examine  the  channel  performance 
for  off-design  conditions. 


3.3  INITIAL  PROFILES 

The  numerical  computation  is  initiated  upstream  of  the  powered 
section  of  the  channel  in  all  cases.  Thus  at  the  initial  station,  MHD 
effects  are  negligible,  and  the  overall  gas  profiles  tend  to  be  similar. 
Given  arbitrary  initial  velocity  and  temperature  profiles,  which  vary 
continuously  between  the  correct  boundary  conditions,  the  profiles 
approach  the  correct  similar  solutions  within  the  first  five  or  six 
x-stations.  Thus  it  is  quite  adequate,  for  the  overall  gas,  to  simply 
impose  parabolic  initial  profiles  and  allow  the  numerical  technique  to 
march  downstream  for  several  integrations  before  entering  the  influ¬ 
ence  of  the  electric  and  magnetic  fields.  The  electron  temperature 
profile  is  equally  simple  to  initiate,  since  the  electron  gas  is  in  ther¬ 
mal  equilibrium  with  the  heavy  species  until  it  enters  the  electric 
field. 

The  initiation  of  the  electron  concentration  profile  is  not  so  simple. 

It  was  found  that  if  the  Ce  profile  is  initiated  in  Saha  equilibrium,  it 
takes  approximately  two  hundred  x-stations,  or  five  feet,  to  asymptoti¬ 
cally  approach  the  correct  Ce  profile.  This  is  because  of  the  relatively 
low  mobility  of  the  ions  in  the  ambipolar  diffusion  process  by  which  the 
Ce  profiles  are  filled  in.  With  this  in  mind,  the  initial  profiles  for  the 
powered  accelerator  section  were  computed  by  starting  the  computer 
program  three  hundred  x-stations,  or  ten  feet,  upstream  of  the  powered 
section.  Three  runs  were  made,  corresponding  to  a  recombination  rate 
of  zero  (frozen  flow),  a  physically  realistic  recombination  rate  as  given 
by  Eq.  (16),  and  a  near -equilibrium  recombination  rate  given  by  Eq.  (16) 
multiplied  times  10^.  The  resulting  initial  Ce  profiles  are  plotted  in 
Fig.  3.  Note  that  on  the  linear  scale  in  Fig.  3,  the  finite -rate  and  the 
near-equilibrium  profiles  appear  identical,  thus  the  boundary  layer  is 
essentially  in  ionizational  equilibrium  for  the  non-MHD  case. 
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The  initial  velocity,  gas  temperature,  and  electron  temperature 
profiles  are  discussed  in  Section  4.1, 


SECTION  IV 
RESULTS 


Computer  calculations  are  discussed  in  detail  in  this  section  for 
the  actual  physical  case,  which  involves  the  simultaneous  presence  of 
electron  thermal  nonequilibrium  and  ionizational  nonequilibrium.  The 
effects  of  varying  the  energy -loss  factor  (6eff)  are  presented,  with  a 
physically  realistic  recombination  rate  coefficient  (ofrec)-  Then  with 
a  physically  realistic  5eff'  the  effects  of  varying  arec  are  explored. 

Nonequilibrium  calculations  are  also  presented  for  the  locally  Hall- 
neutralized  case,  and  compared  with  the  results  of  the  core  flow  Hall- 
neutralized  calculations. 


4.1  SIMULATANEOUS  THERMAL  AND  IONIZATIONAL  NONEQUILIBRIUM 

In  order  to  obtain  a  physically  meaningful  standard  of  comparison, 
the  first  channel  flow  solution  was  run  with  6eff  =  7.  0  x  10 "4  and  orec 
as  given  by  Eq.  (16),  which  were  felt  to  be  the  most  realistic  values 
currently  available.  The  results  of  this  computation  are  presented  in 
Figs.  4  and  5.  Figure  4a  shows  the  geometry  of  the  accelerator  section, 
and  the  configuration  of  the  electric  and  magnetic  fields.  The  nature  of 
the  core  flow  was  primarily  determined  by  fixing  the  magnetic  field  (B), 
the  transverse  current  density  (jy^),  and  the  static  temperature  (T,,,), 

constant  in  the  core  flow  throughout  most  of  the  channel,  with  B  =  2.0 
webers/m2,  jy^  =  15  amps/m2,  and -T^  =  6300°R.  The  pressure  varia¬ 
tion  and  the  area  variation  were  determined  by  Eqs.  (69)  and  (65), 
respectively,  rather  than  being  prescribed  a  priori.  At  the  entrance  of 
the  channel,  there  were  no  fields  present,  with  B  being  increased 
linearly  from  zero  at  x/j?0  =  0  to  its  full  value  at  x/i o  =  10.  Ey  was 
calculated  from  Ohm's  law,  assuming  a  linear  buildup  of  jya>  to  its  full 
value  at  x/j2q  =  10.  Throughout  the  channel,  Ex  was  calculated  by 

,  which  kept  the  core  flow  Hall -neutralized  (jx  =  0). 

00 

) 

The  variations  of  core  flow  pressure,  temperature,  and  velocity  down 
the  channel  are  plotted  in  Fig.  4b. 
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In  Fig.  4  c,  several  of  the  parameters  which  characterize  the 
boundary  layer  are  plotted  down  the  channel.  It  should  be  emphasized 
at  this  point  that  no  attempt  to  satisfy  Maxwell's  equations  has  been 
made.  Therefore,  the  solutions  presented  in  this  investigation  are  not 
necessarily  valid  in  the  region  where  0  <  x/j20  <  10.  Thus  the  rather 
abrupt  variation  of  Cf  in  the  entrance  should  be  disregarded.  The  6 
which  is  plotted  in  Fig.  4c  is  the  boundary  layer  thickness  at  the  point 
where  u/u,,,  =  0.  95.  Note  the  slight  decrease  in  6  and  6*  near  the 
channel  exit.  The  reason  for  this  is  that  boundary  layer  shorting 
developes  near  the  end  of  the  channel,  and  the  resulting  hump  in  the 
velocity  profile  results  in  a  slight  loss  of  mass  in  the  velocity  boundary 
layer.  Local  heat  transfer  rates  are  plotted  in  Fig.  4d.  The  increase 
in  overall  heat  transfer  (q)  near  the  entrance  should  be  ignored,  but 
the  increase  near  the  exit  is  another  significant  result  of  the  B-wall 
shorting  effect.  Note  that  the  heat  transfer  rates  for  the  ion  (q-j)  and 
electron  (qe)  gases,  while  much  smaller,  are  quite  sensitive  to  the 
shorting  phenomena.  It  is  most  significant  that  q0  and  q^  are  several 

orders  of  magnitude  smaller  than  q.  Considerably  higher  operating 
temperatures  would  be  necessary  for  qe  and  q^  to  become  important. 

The  energy  of  recombination  is  arbitrarily  assigned  to  the  ions;  this 
is  why  q^  is  about  an  order  of  magnitude  higher  than  qe. 

Figure  5  illustrates  the  development  of  the  boundary  layer  profiles 
as  they  proceed  down  the  channel.  The  initial  profiles  are  represented 
by  the  first  curve  (x/j2q  =  0)  shown  in  each  figure.  The  most  important 
effect  which  is  illustrated  is  the  presence  of  B-wall  shorting  in  about 
the  last  third  of  the  channel.  This  is  most  evident  in  the  current  density 
profiles  (Figs.  5c  and  d)  which  add  vectorially  to  a  maximum  current 
density  of  122  amps/cm2  at  z/6  =  0.34  and  x/£Q  =  120.  Te  and  Ce 
(Figs.  5a  and  b)  also  display  local  maxima  near  the  point  of  maximum 
current  density.  Thus  Fig.  5  provides  a  picture  of  how  B-wall  short¬ 
ing  tends  to  amplify  very  quickly  once  it  starts  with  j  increasing  Te 
and  Ce  through  joule  heating,  which  in  turn  leads  to  higher  cr  and  hence 
increased  j.  The  influence  of  the  shorting  of  transverse  current  may 
be  seen  in  Fig.  5e,  which  shows  the  resulting  velocity  hump  at  the 
accelerator  exit. 

Although  it  is  not  evident  from  Figs.  4  and  5,  the  dominant  terms 
in  the  electron  energy  equation  for  this  case  are  the  joule  heating  terms 
and  the  collision  term,  with  the  conduction  and  diffusion  terms  becom¬ 
ing  important  only  in  the  region  where  z/6  <0.  05.  One  of  the  reasons 
why  conduction  and  diffusion  are  of  little  importance  is  that  the  electro¬ 
static  sheath  tends  to  insulate  the  electron  gas  from  the  cold,  fully 
catalytic  wall.  To  illustrate  this,  at  x/i0  =  120,  the  wall  temperature 
is  530.  °R,  while  the  electron  temperature  at  the  outer  edge  of  the 
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sheath,  TeQ,  is  1940.  °R.  It  is  also  interesting  to  note  that  at  this  point 
C  =  4.  84  x  10“14,  as  compared  with  a  Saha  equilibrium  concentration 

GO 

at  that  point  of  2.  81  x  10_1®.  While  two  orders  of  magnitude  might 
seem  like  a  great  deal  of  ionizational  nonequilibrium,  this  may  be  put 
in  perspective  by  considering  that  it  is  less  than  0.  01  percent  of  the 
core  flow  electron  concentration. 

In  Section  2.  3,  it  was  pointed  out  that  the  collisionless  sheath 
model  is  invalid  if  CeQ  drops  below  lO-1^  Since  CeQ  =  4.  84  x  10“ 14, 

the  sheath  is  in  fact  collision  dominated  from  a  thermal  standpoint, 
although  still  frozen  for  the  purpose  of  recombination.  In  order  to 
check  the  implications  of  this  violation  of  the  collisionless  sheath 
assumption, 1  the  computer  program  was  run  with  TeQ  and  Ceo  set 

equal  to  their  equilibrium  wall  values,  530°R  and  3.  29  x  10“39J  respec¬ 
tively.  Surprisingly  enough,  the  conduction  and  diffusion  terms  were 
still  not  strong  enough  to  cause  any  significant  change  in  the  results  in 
Figs.  4  and  5,  outside  a  z/6  value  of  about  0.2.  Thus,  although  the 
simple  collisionless  sheath  model  leaves  much  to  be  desired,  it  is 
adequate  for  the  calculations  presented  herein. 


4.2  ELECTRON  THERMAL  NONEQUILIBRIUM 

To  isolate  the  effects  of  electron  thermal  nonequilibrium,  the  com¬ 
puter  program  was  run  for  the  same  case  as  in  Section  4.  1,  except  with 
different  values  of  energy -loss  factor.  For  this  discussion,  it  is  useful 
to  use  the  species  energy -loss  factor,  6S,  which  is  related  to  the  effec¬ 
tive  energy-loss  factor,  6eff,  by  Eq.  (33).  The  physically  realistic 
energy-loss  factor  is  6eff  «  7.0  x  10”4,  which  corresponds  to  »  18. 
The  lowest  value  of  which  is  reported  in  the  literature  is  given  as 
5N2  584  10* ,  in  Sutton  and  Sherman  (Ref.  3).  This  lower  energy -loss 

factor  leads  to  much  higher  electron  temperatures,  and  so  a  case  was 
computed  with  this  value.  For  the  equilibrium  extreme,  was 

multiplied  by  a  factor  of  104  to  yield  6^2  «  1.8  x  10^. 

Figure  6  compares  the  boundary  layer  profiles  for  6N2  =  10,  18, 
and  1.  8  x  10^,  at  x/i0  =  82.  0.  Smaller  values  of  fijjg  tend  to  severely 

aggravate  the  B-wall  shorting  problem,  while  no  shorting  at  all  exists 
for  the  case  of  electron  thermal  equilibrium  (6N2  =  1-8  x  10^).  Keep 

in  mind  that  these  three  cases  were  all  run  with  the  same  recombination 
rate  coefficient  given  by  Eq.  (16). 
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Figure  7  shows  the  effect  of  varying  the  energy-loss  factor  on  two 
of  the  most  important  overall  boundary  layer  parameters,  6*  and  q.  In 
both  Figs.  7a  and  b,  the  6^2  =  18  line  is  the  same  as  in  Figs.  4c  and  d. 

The  plots  for  6N2  =  10  31,6  terminated  beyond  x/iQ  »  85  because  the  elec 
tron  temperature  exceeded  the  range  of  validity  of  some  of  the  assump¬ 
tions  used  in  evaluating  the  transport  properties. 


4.3  IONIZATIONAL  NONEQUILIBRIUM 

With  a  realistic  energy -loss  factor,  the  computer  calculation  of 
Section  4. 1  was  repeated  twice  more.  This  time  with  arec  =  0,  to 
duplicate  frozen  recombination,  and  also  with  the  arrec  given  in  Eq.  (16) 
multiplied  times  10 5 ,  to  duplicate  near -equilibrium  recombination.  The 
resulting  profiles  are  displayed  in  Fig.  8. 

Note  that  Ce  does  remain  very  near  to  Saha  equilibrium  for  all 
values  of  z/6  above  0.  3.  Recombination  is  essentially  frozen  then 
from  about  z/6  =  0.  25,  on  in. to  the  wall.  The  frozen  recombination 
profiles  display  no  shorting  because  primarily  a  is  no  longer  coupled 
to  Te  through  Ce.  The  streamwise  variation  of  the  boundary  layer 
parameters  is  not  plotted  since  the  physical  case  is  so  near  the  Saha 
equilibrium  case. 


4.4  LOCAL  HALL-NEUTRALIZATION 

It  is  not  uncommon  practice  to  neglect  the  effects  of  Hall  currents 
when  performing  an  analysis  of  nonequilibrium  effects.  This  corre¬ 
sponds  to  local  Hall -neutralization,  where  Ex  at  every  point  in  the  flow 
is  assumed  to  have  the  value 

j  B 

Ex  =  (71) 

nee 

such  that  no  Hall  current  flows  anywhere.  Then  Ohm's  law  becomes 

jy  -  cr(Ey  -  uB)  (72) 

Such  a  procedure  tends  to  overestimate  the  importance  of  electron 
thermal  nonequilibrium,  as  is  shown  in  Figs.  9  and  10.  Here  the  com¬ 
puter  results  of  Section  4. 1  are  compared  with  a  similar  computation 
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using  Eqs,  (71)  and  (72).  The  reason  for  this  behavior  may  be  seen 
from  the  simplified  electron  energy  equation 


T)  =  Evj 


xJex 


E 


y3ey 


-  n  euE„  + 
e  x 


en 


(73) 


■which  results  from  Eq.  (5)  when  all  but  the  collision  and  joule  dissipa¬ 
tion  terms  are  neglected.  Consider  the  right  hand  side  of  Eq.  (73)  for 
the  core  flow  Hall -neutralized  case.  Here  Eyje  is  the  largest  term,  it 


is  positive  because  both  Ey  and  je  are  positive.  Ex  and  u  are  positive, 

while  jex  is  negative.  Thus,  all  three  of  the  other  terms,  the  Hall 
terms,  on  the  right  hand  side  are  negative.  So  while  the  transverse 
current  j6y  is  putting  energy  into  the  electron  gas,  the  Hall  terms  are 


taking  energy  out.  So,  in  effect,  there  is  some  Hall  generation  going 
on  in  the  boundary  layers  of  this  Faraday  accelerator.  If  the  accelera¬ 
tor  is  run  with  local  Hall -neutralization,  then  two  of  the  negative  Hall 
terms  drop  out,  and  Eyje  becomes  much  larger.  Thus  increasing  the 


joule  dissipation  into  the  electron  gas  to  an  unrealistically  high  value. 


Joule  dissipation  profiles  are  plotted  for  both  types  of  Hall- 
neutralization  in  Fig.  9c.  This  form  of  the  joule  dissipation  term  may 
be  developed  from  Equation  (73)  and  using  Ohm's  law. 
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The  effects  of  local  Hall-neutralization  on  the  overall  boundary  layer 
parameters,  6*  and  q,  are  shown  in  Fig.  10.  Note  that  the  severe 
boundary  layer  shorting  brought  on  by  local  Hall -neutralization  causes 
considerable  overestimation  of  the  heat  transfer.  Again,  calculations 
for  the  locally  Hall-neutralized  case  were  terminated  beyond  x/iQ  =  35.  2 
because  of  the  high  electron  temperatures  involved. 


SECTION  V 
CONCLUSIONS 


For  the  class  of  nitrogen  accelerators  under  study  here,  nonequilib¬ 
rium  effects  are  of  secondary  importance  in  laminar  boundary  layers 
to  an  x/ SL0  of  about  ninety,  but  beyond  that  point  they  are  a  direct  cause 
of  B-wall  shorting.  Once  shorting  begins,  electron  thermal  nonequilib¬ 
rium  causes  it  to  amplify  very  quickly,  causing  significant  increases  in 
heat  transfer.  If  the  channel  is  designed  for  equilibrium  boundary 
layers,  then  the  B-wall  shorting  is  sufficient  to  cause  thermal  choking 
near  the  exit.  Although  there  is  considerable  ionizational  nonequilib¬ 
rium  in  the  colder  portions  of  the  boundary  layer,  it  is  insufficient  in 
the  hotter  areas  (where  appreciable  current  flows)  to  affect  the  solutions 
noticeably.  The  relatively  simple  model  of  Hale  and  Ker  re  brock  (Ref.  1), 
which  uses  only  the  collision  and  joule  heating  terms  in  the  electron 
energy  equation  and  which  assumed  Saha  equilibrium  at  the  electron 
temperature,  would  have  been  quite  adequate  to  describe  the  nonequilib¬ 
rium  effects  for  this  application.  Hall  effects  can  be  very  important  to 
the  overall  operation  of  a  core  flow  Hall -neutralized  device  and  should 
be  neglected  only  if  there  is  qualitative  justification  for  it. 

This  study  also  points  out  the  need  for  improved  understanding  of  the 
basic  collision  phenomena  involved  in  plasmas  of  engineering  interest. 

For  example.  Figs.  6  and  7  demonstrate  that  the  boundary  layer  behavior 
is  very  sensitive  to  the  value  of  the  energy-loss  factor.  While  reliable 
energy-loss  data  are  abundant  for  low  temperature  nitrogen  and  relatively 
high  energy  electrons,  it  is  very  scarce  for  higher  temperature  gases 
involving  lesser  degrees  of  thermal  nonequilibrium.  Thus  it  should  be 
emphasized  that  these  boundary  layer  calculations  are  only  as  good  as 
the  cross  sections  and  energy-loss  factors  involved,  and  any  future  work 
undertaken  should  be  based  on  improved  evaluations  of  these  quantities 
as  they  become  available. 
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a.  Geometry  and  Field  Configuration 

Fig.  4  Streamwise  Variation  of  Channel  Flow  Parameters  using  Physically  Realistic 
Recombination  Rate  and  Energy- Loss  Factor 
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b.  Core  Flow  Properties 
Fig.  4  Continued 
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c.  Boundary  Layer  Parameters 
Fig.  4  Continued 
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T  and  Te  (°R  x  io3) 

a.  Overall  and  Electron  Gas  Temperature 

Fig.  5  Development  of  Boundary  Layer  Profiles  using  Phsically  Realistic  Recombination  Rate  and  Energy- Loss  Factor 


d.  Hall  Current  Density 
Fig.  5  Continued 
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e.  Velocity 
Fig.  5  Concluded 
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c.  Transverse  Current 
Fig.  6  Concluded 
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b.  Local  Heat  Flux 
Fig.  7  Concluded 
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c.  Transverse  Current  Density 
Fig.  8  Concluded 
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b.  Trasverse  Current  Density 
Fig.  9  Continued 
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a.  Displacement  Thickness 

Fig.  10  Streamwise  Variation  of  Boundary  Layer  Parameters  for  Local  and  Core  Flow  Hall-Neutralization 
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